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ABSTRACT 

The auatomical structure of the braiu cau be observed via 
uou-iuvasive techuiques such as diffusiou imagiug. How¬ 
ever, these are imperfect because they miss couuectious that 
are actually kuowu to exist, especially loug rauge iuter- 
hemispheric oues. lu this paper we formulate the iuverse 
problem of iuferriug the structural couuectivity of braiu uet- 
works from experimeutally observed fuuctioual couuectivity 
via fuuctioual Maguetic Resouauce Imagiug (fMRI), by for- 
mulatiug it as a couvex optimizatiou problem. We show that 
structural couuectivity cau be modeled as au optimal sparse 
represeutatiou derived from the much deuser fuuctioual cou¬ 
uectivity iu the humau braiu. Usiug ouly the fuuctioual 
couuectivity data as iuput, we preseut (a) au optimizatiou 
problem that models coustraiuts based ou kuowu physiologi¬ 
cal observatious, aud (b) au ADMM algorithm for solviug it. 
The algorithm uot ouly recovers the kuowu structural cou¬ 
uectivity of the braiu, but is also able to robustly predict 
the loug rauge iuterhemispheric couuectious missed by DSI 
or DTI, iucludiug a very good match with experimeutally 
observed quautitative distributious of the weights/streugth 
of auatomical couuectious. We demoustrate results ou both 
syuthetic model data aud a hue-scale 998 uode cortical dataset, 
aud discuss applicatious to other complex uetwork domaius 
where retrieviug effective structure from fuuctioual sigua- 
tures are importaut. 

General Terms 

Theory, Algorithms, fMRI, structural, fuuctioual, couuec¬ 
tivity, uet works 
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Networks or graphs are used to model the structure aud 
dyuamics of mauy real world systems, such as the iuter- 
uet aud the world wide web, scieutihc collaboratious, iu- 
frastructure systems such as road, rail, power or trausport 
uetworks, the braiu, social uetworks of iudividuals, kuowl- 
edge aud desigu structures iu eugiueeriug, iuformatiou aud 
kuowledge dyuamics iu hrms or orgauizatious 
[^. Usually, au implicit assumptiou is that the structure of 
the uetwork is observable, aud forward models or iufereuce 
of dyuamics cau be based ou the topological orgauizatiou of 
the system, for example [^. 

However, mauy of these systems, uotably the braiu, have 
structures that are difficult to map but cau be observed 
via the activity they support. Further, kuowledge of fuuc¬ 
tioual dyuamics cau also provide deeper iuformatiou ou ac¬ 
tive structure (as a subset of the total structure) aud its 
temporal variatiou iu systems. Therefore, we cau dehue the 
inverse problem of deriving strueture from the funetional sig¬ 
nature of dynamieal networks. The auatomical structure of 
braiu uetworks cau ouly be observed via techuiques such as 
diffusiou imagiug, [e.g. Diffusiou Spectrum Imagiug (DSI), 
Diffusiou Teusor Imagiug (DTI)], which are imperfect be¬ 
cause they miss edges that are kuowu to actually exist iu 
the system, especially loug rauge oues betweeu the two hemi¬ 
spheres [^[^[^[^1^. Further, uot all of the auatomical 
couuectious are active at all times; some are dormaut, lead- 
iug to the uotiou of effeetive eonneetivity. It is possible to ob¬ 
serve fuuctioual couuectivity iu these systems, measured iu 
terms of correlatious of uetwork activity via techuiques such 
as fuuctioual Maguetic Resouauce Imagiug (fMRI), EEG, 
or MEG. Methods that cau iufer the auatomical, structural, 
or effective couuectivities, startiug from the fuuctioual sig- 
uature will therefore provide a siguificaut step forward iu 
uuderstaudiug the uormal aud/or diseased structural aud 
fuuctioual states of the braiu by providiug iusight iuto the 
liuks betweeu structure aud fuuctiou. 

Solviug either the forward problem (iuferriug dyuamics 
from structure) or the iuverse problem (iuferriug structure 
from dyuamics) exactly is curreutly beyoud the capacity of 
auy discipliue. However, receut efforts iu the braiu uetworks 
domaiu show that both the forward aud iuverse problems 
are extremely topical iu several discipliues 
across computatioual ueuroscieuce, data miuiug aud physics. 
Iu this paper, we pose the iuverse problem of miuiug the 
struetural eonneetivity of braiu uetworks from a siguature 
of its funetional eonneetivity. Differeut approaches iu dif- 
fereut discipliues model the problem iu several ways. Eor 
example, computatioual ueuroscieuce relies heavily ou ex- 


periments followed by the use of established graph theoretic 
analysis , whereas physics based approaches rely on the¬ 
oretical modeling of physiological phenomena, such as the 
expected behavior of eigenvalues of systems in stable states 
or near critical regimes [^[^|^. In data mining, recent ap¬ 
proaches have focussed on deriving functional connectivity 
signatures from functional connectivity data, using sparse 
representations |18[ , but we are not aware of any work 

on deriving structural connectivity from functional connec¬ 
tivity using sparse representations. 

This paper addresses some significant gaps that remain 
unexplored in previous (very recent) attempts on model¬ 
ing the inverse problem [21[ [^. Firstly, we focus on the 
sparse structure of structural connectivity as opposed to the 
dense structure of functional connectivity: an optimal solu¬ 
tion will be sparse, since experimental measurements show 
that the fMRI signature is always much denser than the 
structural/anatomical network. Secondly, while it is reason¬ 
able to assume that functional connectivity ultimately de¬ 
rives from structure and so the two would share topological 
properties, we focus on the issue of the actual quantitative or 
numerical distributions of the connection strengths in both 
cases. As experimental data shows, these are very different 
for the structural and functional data, with the structural 
DSI or DTI data being non-negative and the fMRI data 
containing both positive and negative correlations. Addi¬ 
tionally, the numerical ranges and the distribution of the 
data are very different for both cases. 

In this paper, we propose a model for deriving a sparse 
weighted representation of the structural connectivity of the 
brain from its dense (full matrix) weighted functional con¬ 
nectivity. Our main contributions are as follows: 

1. Formulation: A convex optimization problem dehni- 
tion of deriving the structural connectivity of the brain 
from its fMRI connectivity. 

2. Algorithm: An ADMM based algorithm for solving 
the optimization problem. 

3. Application: Applying the approach to real experi¬ 
mental DSI and fMRI data. 

4. Prediction and validation: To what extent does the 
inferred structural connectivity confer with the known 
structural connectivity in the cortex and can it pro¬ 
vide any new information about anatomical connec¬ 
tivity? We compare our results with (a) the known 
short range anatomical connectivity of the cortex, and 
(b) the known (but missing from experimental data) 
long range inter-hemispheric connectivity of the cortex 
(prediction of missing links in the structural network). 

2. BACKGROUND 

Before formalizing the problem, we present some obser¬ 
vations that will form the motivation for our optimization 
formulation, and also establish the basis for why posing the 
problem within an optimization framework can have signif¬ 
icant advantages: 

1. Experiments have shown that the structural and func¬ 
tional connectivity of the brain have high positive cor¬ 
relations with each other [^. In other words, if two 
areas of the brain are strongly connected anatomically. 



Figure 1: ^lExperimentall^ measured anaffbmical and 
functional connectivity in the human cortex [1Q| . (a) 
Structural connectivity measured via DSI. (b) Rest- 
ing state functional connectivity measured via fMRI. 
(c) Linear correlation map between structural and 
functional connectivity matrices. 


it is likely that they will share high functional corre¬ 
lations too. Figures a) and (b) show the anatomical 
and resting state functional connectivity data from . 
Each row or column of the matrix shows a region of 
interest or ROI of the cortex, and its entries represent 
structural or functional connectivity with all the other 
ROIs. The matrices represent both the left and right 
hemispheres, with the top half of the rows/columns 
representing one hemishphere and the bottom half the 
other. Figure[^c) shows the linear correlation between 
the rows/columns of the structural and functional con¬ 
nectivity matrices. 


2. Experiments show that the functional connectivity is 
based upon direct as well as indirect connections or the 
presence of both direct and indirect structural paths 
between brain regions 21 Thus, signatures 

of functional connectivity or correlations of activity as 
mapped via fMRI are significantly denser than signa¬ 
tures of structural connectivity as mapped via tech¬ 
niques such as DSI or DTI [^, since the functional 
connectivity is postulated to derive from both direct 
and indirect connectivity between brain regions. As 
can be seen from Figure [^b), the functional signa¬ 
ture is much denser than the structural one in Fig¬ 
ure Ha). Thus, it appears from experimental data that 
structural connectivity can be naturally expressed by 
a sparser representation while the functional one is a 
denser one. 


3. Current non-invasive techniques of mapping structural 
brain connectivity miss the long range connections, es¬ 
pecially inter-hemispheric ones in the brain [19[ 21 ^ 
10 . This is visible in Figs[^a) and (b): while the main 


diagonal is similarly structured in both, the two off- 
diagonals, representing inter-hemispheric connectivity 
between the two hemispheres, are very weak in Fig ^a) 
and strong in Fig[^b), showing that the functional sig¬ 
nature captures traces of long range connections via 
strong correlations that the structural data from DSI 


The above observations lead to the following informal char¬ 
acterization of the problem. Given a dense functional con¬ 
nectivity matrix, can we derive an optimal sparse repre¬ 
sentation of structural connectivity, by modeling physio¬ 
logically and experimentally known constraints on struc¬ 
tural and functional connectivity? Further, (a) how much 
of experimentally mapped and known structural connectiv¬ 
ity is the sparse representation able to predict, and (b) can 

























the sparse representation also predict structural connectivity 
that is known to exist anatomically, but is missing from data 
that comes from the current state of the art non-invasive 
imaging techniques? 

3. PROBLEM DEFINITION 

Formally, we consider two graphs Gs = (V, Es) and G/ = 
(y,Ef), with ROIs of the cortex represented as N nodes 
in the set V, and the edge sets Eg and Ef representing 
structural and functional connectivity edges, respectively, 
between two ROIs or nodes. While structural connectivity 
is a measure of the actual anatomical connectivity between 
two ROIs or nodes, functional connectivity is a measure of 
correlation between activity on two nodes as the brain is in 
resting state or performs a specific task. Thus, note that 
while the presence of an edge in the set Eg is a measure 
of actual physical connectivity, an edge in the Ef signifies 
the level of correlation between activity based on the BOLD 
signature [^ . 

However, based on the observations in the background 
and previous research , we 

assume that the structural and functional connectivity of a 
node to other nodes is positively correlated, as discussed in 
the Background section. More formally, we represent the ex¬ 
perimentally measured structure and functional connectivity 
networks Gg and Gf hy weighted adjacency matrices S and 
E respectively, where the entries Sii and Eij represent the 
structural and functional connection weights between brain 
regions i and j, respectively. Thus, each column (or row) Si 
or fi represents the structural or functional connectivity of 
node i to all the other nodes in the network. 

Our main task can be stated as inferring a sparse rep¬ 
resentation of T, as a representation of the structural con¬ 
nectivity, by modeling constraints that derive from known 
experimental and physiological observations. The simplest 
way of stating this is 

FX = F, ( 1 ) 

where X represents our sparse representation matrix, and we 
model it as a transformation that takes F to itself. However, 
to extract the required pattern, this form can be further sim¬ 
plified to our advantage if we work with a lower dimensional 
representation of E rather than the full N x N version. The 
reasoning for this is as follows. Row i of E shows the connec¬ 
tivity of node i with all the other connections. Because of the 
physiologically known modularity structure of the brain (e.g. 
visual cortex, auditory cortex, etc.), we know that several 
other nodes will typically have connectivities very similar 
to node i. Thus, we are looking for a representation where 
each node can be represented as a point in /c—dimensional 
space, with k « N, such that if two nodes share similar 
connectivities to other nodes in T or S', then they must lie 
very close to each other in the /c—dimensional space. This 
can be achieved in several ways, but here we use the sim¬ 
plest possible representation: the spectral representation of 
T, or the co-ordinates given by the first k eigenvectors of E; 
E = VDV^, Ek = VkDkV^ , where 14 which iskxN matrix 
, now represents the positions of N nodes in /c—dimensional 
space R^. When the transformation X is applied to T, the 
components of sparse X will be like weights on the columns 
of Vk , “picking out” the most relevant connections of node i 
to all other nodes. Thus we have 

wx = w (2) 


This is also the simplest possible representation for which 
the known physiological constraints described in the previ¬ 
ous section can be captured: (a) high correlations in the 
fMRI data must imply high probability of structural con¬ 
nectivity, but (b) the structural connectivity must be sparse, 
since only direct connections must be inferred, and (c) both 
the short range and the long range structural connectivity 
must be inferred using the functional connectivity. 

One complication is introduced by the presence of negative 
entries in the experimental fMRI data. The meaning of neg¬ 
ative correlations in fMRI data is under open debate and has 
been attributed to (a) pre- and post-processing of the fMRI 
data (especially removals related to other physiological func¬ 
tions such as cardiovascular and respiratory functions from 
the raw fMRI data, and removals of global modes when tem¬ 
poral averaging is performed to reveal spatial correlations as 
demonstrated in Fig.I^b), (b) BOLD activity and haemody¬ 
namic effects, and (cjexcitatory versus inhibitory effects as 
correlated with positive or negative BOLD response. Since 
the exact meaning of the negative correlations in the fMRI 
matrix and its physiological significance is unclear and var¬ 
ied in the research literature, and all structural connectivity 
data as measured by DSI or DTI techniques is positive, we 
would like to separate out the positive and negative parts, 
as 

W(Xp + Xn) = W, (3) 

where Xp is sparse and positive, and Xn is small and neg¬ 
ative. Separating out the positive and negative also makes 
our results easier to interpret. Surprisingly, we discover in 
our results that both the positive and negative solutions echo 
back a similar “picture” of structural connectivity, but Xp 
is much closer to predicting the real structural connectiv¬ 
ity. The robustness of this result is particularly significant, 
because unlike previous efforts [^, we make no use of the 
structural connectivity matrix to predict the same in our 
optimization formulation; all our computations are done on 
the functional connectivity matrix and what we derive is 
then shown to be very close to the experimentally measured 
structural connectivity. 

Thus, the final optimization problem can be stated as 
min||Xp||i + ^||X„||2 

sub to Vk{Xp + Xn) = Vk 

diag{Xp) = 0 ( 4 ) 

diag{Xn) = 0 

Xp > 0 

> 0. 

Note that, similar to [^, we want to enforce the diagonals 
of the solution variables to be 0 , so as to avoid the triv¬ 
ial solution of each node expressing itself as its own linear 
combination and none of the others. 


4. OPTIMIZATION ALGORITHM 

We present an ADMM based solution framework for solv¬ 
ing the above problem. Introducing auxilliary variables A 
and B for the optimization variables Xp and Xn, and indi¬ 
cator functions /i+(Xp) and / 2 +(—Xn) for the non-negative 


and negativity constraints, we have: 

L[A,B,X^,Xn,Ai,A2] = 

II^pIIi + - iA + B)V^\\l 

+tr[Af - (Xp - diag{Xp))] 

+tr[Al{B - {Xn - diag{Xn))] (5) 
+ ^\\A-{Xp-diag{Xp))\\l 
+ Pl\\B - {Xn - diag{Xn))\\l 
+h+{Xp) + h+i-Xn), 


4.4 Updating Xn 

Similar to the derivation for A and B, we can differentiate 
with respect to X„ as 

dL _ d 
dXn ~ dXn 

+ ^\\B - Xn + diag{Xn)\\l (14) 

d-tr [A 2 {B — Xn + dicig{Xn'}'}] , 

set the derivative to 0 and apply indicator function / 2 +(—Xn) 
to get 



where /i+(Xp) = 0 when Xp > 0 and oo otherwise, and 
/ 2 +(—Xn) = 0 when Xn < 0 and oo otherwise. 

Now we minimize L (w.r.t. its arguments) using a stan¬ 
dard iterative ADMM process where we differentiate with 
respect to one variable while keeping the others fixed and fi¬ 
nally updating the Lagrange multipliers Ai and A 2 . We rep¬ 
resent with A*, B*, X*, X*, Ai, A 2 , the updated variables 
A, B, Xp, Xn, Ai, A 2 . The full algorithm is presented in Al¬ 
gorithm 

4.1 Updating A 

A can be updated as computing 


dA 


dA 


^\\Vk-{A + B)Vk\\l 


+ ^||j4 — Xp + diag{Xp)\'^ 


-\-tr [A^ {A — Xp + diap(Xp))] , 


( 6 ) 


and setting it to 0. We then get 

V = 14+Pi/)“'[A*if’(t4-B14) + piXp-Ai]. (7) 


4.2 Updating B 

The form for B is exactly similar to A. Thus, 

B* = {\tV^ 14+P2/)“'[Atlf {Vk-AVk)^P2Xn-A2]. (8) 


4.3 Updating Xp 

Since the norm of Xp is to be minimized, along with 
the square terms, Xp can be updated in closed form using a 
soft-thresholding operator as follows: 

X;=Y-diag{Y), (9) 

where 


Y = Tj_(a+—), (10) 

PI V pi J 



with 


Tp{v) = (|i;| - ri)+sgn{v). 
Now applying the indicator function /i+(Xp), 


X; = Y -diag{Y), 



where Y takes only positive values due to /i+(Xp). 


( 12 ) 


(13) 


where I is the identity matrix and Xn takes on only negative 
values because of / 2 +(—Xn). 


4.5 Updating Ai and A 2 

Finally, we update the Lagrange multipliers as 


At = Ai +pi(V -x;), 

(16) 

A^ = A2+P2(B*-X:). 

(17) 


Algorithm 1 ADMM Algorithm 


Input: F,k,Xt,pi,p 2 
Output: Xp^Xn 

Initialize Xp, Xn, A, B, Ai, A 2 

1: VDV^ 

2: Use 14 from VkDkVj^ (first k eigenvectors) 

3: repeat 


4: 

5: 

6 : 

7: 


9: 

10 : 

11 : 


A ^ (A*lf 14+pi/)-MAtlf (14-P14)+piXp 
B ^ (Atlf 14+p2/)-MAtlf (Vfe-4Vfe)+p2X„ 


Y 

V 


A+ Ai- AiiT 

PI PI 

p^Y- diag{Y) 




BFA 2 


Xn ^ Z — diag{Z) 

Ai^Ai+pi(A*-Xp*) 

A2^ A^ = A2 + P2(B*-X;) 


12: until convergence 


13: Xn 


14: Xn 




Ai] 

A 2 ] 


5. RESULTS 

5.1 Synthetic baseline model -1 

As discussed in the Background section, the real brain 
data that we will work with has known missing entries, es¬ 
pecially for long range interhemispheric connections. In or¬ 
der to test and validate the approach on a baseline, we have 
generated synthetic models mimicing the known large scale 
architecture of the cortex. Studies have shown that the 
cortex has a hierarchical modular architecture with larger 
modules nesting smaller modules on several levels 
|15] . We have therefore generated hierarchically modular 
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Figure 2: Synthetic baseline model I. (a) adjacency 
matrix for a hierarchical network with brain like ar¬ 
chitecture. (b) derived synthetic functional matrix, 
(c) Xp (d) Xn (e) Xpt with near zero entries removed, 
(f) Xpn ^XpH {Xn = 0). Parts (c), (d), (e), and (f) 
all bring out the connectivity structure of the adja¬ 
cency matrix from the functional matrix. Although 
the hierarchical structure is not particularly visible 
in (e) and (f), it is detected, as confirmed by preci¬ 
sion and recall measurements. 



k 

Figure 3: Precision and recall with respect to 

varying number of dimensions preserved in spec¬ 
tral representation k for synthetic baseline model 
I. Solid lines show precision, dotted lines show re¬ 
call, for Xp, Xpt with near zero entries removed and 
Xpn = Xp n {Xn = 0). Xp has low precision and almost 
perfect recall. Precision climbs significantly for Xpt 
and Xpn for lower recall of upto about 0.80. 


networks with stochastic block model type architectures, fol¬ 
lowing closely the models introduced in ^3 1^. This is the 
structural connectivity matrix S. Figure^ a) shows an ex¬ 
ample hierarchical network adjacency matrix S, with 1024 
nodes, and 16 modules at the finest scale along the main 
diagonal. Other hierarchical levels are formed by progres¬ 
sively reducing the probability of connectivity for 8 mod¬ 
ules, 4 modules, and finally 2 modules representing the left 
and right hemishpheres. Further, a similar architecture is 
repeated for inter-hemispheric connectivities, appearing as 
the off-diagonals in the matrix. 

In order to generate a similar synthetic functional activ¬ 
ity / correlations network, we have used the idea already 
observed in the Background section and Fig.[^ that the ex¬ 
istence of direct as well as indirect paths in the structural 
matrix are highly positively correlated with measured func¬ 
tional activity in the cortex. That is, if two brain regions 
are connected by direct or indirect paths, then it is highly 
likely that their functional activity correlation signature will 
also be high. Shorter the path, more direct the connection, 
higher the probability for positive functional correlations. 
Following previous work done in [21[ |19| and experimental 
observations in [12[ , a simple model of producing a syn¬ 

thetic functional matrix from a given structural matrix is 
defined as: 

F = S + S^ + ... + S^, (18) 

where each power matrix captures the numbers of paths of 
length j in the adjacency matrix. Analytically j can go upto 
infinity, but practically and computationally, a small number 
can be chosen to simulate this. The reasoning behind this 
is twofold: (a) it can be assumed that higher powers (longer 
and longer paths) contribute less and less to the functional 
activity, and (b) the brain operates at marginal stability, 
close to but lower than the critical boundary between stable 
and unstable, thus implying the condition that the principal 
eigenvalue of S will always be less than 1. In such a case, 
the series in Eqn converges, as shown in [^. For more 
analytical detail, please refer to previous work in . 

Figure l^b) shows a synthetic functional connectivity ma¬ 
trix derived from the structural connectivity matrix by sum¬ 
ming upto paths of length 5. This is matrix of functional 
connectivity, denoted by F. 

We apply our algorithm onto the spectral representation 
of F, and perform an analysis on the structural connectiv¬ 
ity extracted when we vary k, the number of dimensions 
preserved in the spectral representation. Figures l^c) and 
(d) show the matrices Xp and Xn, respectively at A; = 16. 
Further, Xp shows a sparse structure with several entries 
extremely close to zero (by at least two or more orders of 
magnitude). We remove these near zero values, producing 
the matrix Xpt, Fig. [^e), that brings out the structure of 
the adjacency matrix. Further, as is seen in Figs [^c) and 
(d), both Xp and Xn bring out the structure in two different 
ways: the positive entries in Xp bring out structure, and the 
zero and near-zero negative entries of Xn also mimic and re¬ 
inforce the same structure, since they imply the absence of 
negative relationships between the nodes. In the case of the 
synthetic baseline we work with here, the negative entries of 
Xn are quite small, since there are no negative entries in F. 
But as we will see in the next section, Xn serves an impor¬ 
tant role when F contains negative entries. Thus, we can 
combine information by extracting Xpn = Xpt fl {Xn = 0), 














which is shown in Fig. [^f). All these interpretations bring 
out similar architectures for the adjacency matrix S, start¬ 
ing from the functional matrix F. Observe that we make no 
use of S to derive Ap, Apt, and Ap^. 

Figure!^ a) shows the precision and recall for Ap, Apt, and 
Apn. Precision is the number of correctly identified connec¬ 
tions in our solution divided by the total number of identified 
connections. Recall is the number of correctly identified con¬ 
nections divided by the total actual number of connections 
that actually exist in the structural connectivity matrix and 
should be ideally identified. Ap has low precision and al¬ 
most perfect recall. Precision climbs significantly for Apt 
and Apn for lower recall of upto about 0.80. This observa¬ 
tion will become important in the next section, because for 
real brain data, precision will drop while recall will remain 
high. The reason for this is that we will predict several miss¬ 
ing links with the algorithm, notably the inter-hemispheric 
connectivity, which will necessarily bring down the preci¬ 
sion. Therefore, it was important to check with a baseline 
model that had all the relevant connectivity been present 
in the structural data, both precision and recall would have 
had acceptably high values. 


5.2 Brain data and data-driven baseline - II 

The brain data that we work with to test our algorithm 
comes from . A primary reason for choosing this dataset 
is that it is the only dataset we could find that maps struc¬ 
tural (DSI) and functional (fMRI) connectivity for a defined 
set of 998 region parcellation map of the human cortex. 
While there are other datasets available that contain either 
only structural connectivity or only functional connectivity, 
it is often hard to put together these data and transfer re¬ 
sults of one dataset to another since the parcellations used 
may be different in every case. These problems are discussed 
in [^. A second reason is that our results and methods 
could be tested against that also use the same dataset 

for the same aim as ours. Thirdly, while there is work on 
task based or disease based functional networks with tem¬ 
poral slices of functional connectivity, for e.g. 25 m, the 
problem that we were working on involves the use of tem¬ 
porally averaged spatial correlations between brain regions, 
such that the structural/physical/spatial connectivity can 
be derived. For this reason, we have used temporally av¬ 
eraged, resting state fMRI data, along with a parallel map 
of DSI based structural connectivity data for the same set 
of 998 regions of the cortex. Figure [^a) and (c) show the 
structural DSI and functional fMRI data, respectively, with 
the histograms of their connections shown in Figs[^b) and 
(d). Howeve, due to the general nature of the way we have 
formulated the problem, our algorithm could easily be ap¬ 
plied to temporal data networks or disease networks, with 
the expectation that the effective structural connectivity or 
abnormal structural connectivity, respectively, could possi¬ 
bly be retrived from functional signatures of task based fMRI 
data or diseased brain fMRI data. 

The first baseline question we now ask is: if structure and 
functional data show high positive correlation, as claimed 
by the literature and discussed in the Background section, 
and is somewhat confirmed by a superficial visual look at 
the two datasets, can we apply a much simpler threshold¬ 
ing algorithm to derive the structural connectivity from the 
functional one? We tested this using Algorithm 2. 

The results are shown in Fig.j^e) and (f). Note that we 



Figure 4: Brain data from [10| and co mpar isons 

with baseline solution and solution from |2l[ . (a) 

Structural connection matrix, (b) Histogram for 
(a), (c) fMRI functional connection matrix, (d) 

Histogram for (b). (e) Simple thresholding based 

baseline solution derived from (c). (f) Histogram 

for (e) showing that simple thresholding may bring 
out topological connections weakly, bu t fa ils on con¬ 
nection strengths, (g) solution from [21| . (h) his¬ 

togram for (g): solution in [21| correctly retrieves in- 
terhemispheric connectivity with largest eigenvalue 
of derived structural connectivity = 0.87 that re¬ 
spects stability criteria, but is not sparse, has nega¬ 
tive entries, and has a different range of connection 
strengths as compared to (b). 






































Algorithm 2 Simple Thresholding Algorithm 

Input: F, S 
Output: F* 

1: z ^ Compute the number of non-zero entries in S 
2: o ^ Order entries of F (+ve and -ve) by magnitude in 
descending order 
3: t ^ o(z) 

4: F* ^ (|F|>t) 

5: Compare properties of S and S* 


cannot compute precision and recall in any reliable way for 
the derived F*, because we have used the number of non¬ 
zero entries in S to decide the threshold t. The results show, 
however, that even though the visual form of F* in Fig.l^e) 
does have some similarity to S in Fig.[^a), the histograms of 
connection strengths look quite different. This happens be¬ 
cause thresholding by absolute value to bring the sparsity of 
F* close to the sparsity of S causes values to intermittently 
disappear from F, since F contains both negative and pos¬ 
itive values. Further, this simple thresholding mechanism 
uses information in F, which our algorithm does not. 

Finally, we also study results from III and [^. The re¬ 
sults in have used the coarser 66 x 66 region definition, 
and are therefore not as fine grained as the 998 x 998 re¬ 
gion definitions used in our work (the same source data is 
available at both resolutions, but using the higher resolu¬ 
tion is more challenging because it introduces much more 
fine structure and constraints into the problem). The re¬ 
sults from are shown in Fig.j^g). The biggest strength 
of this result is that there is a physiologically based model 
that explains both the forward and inverse problems, along 
with predicting the missing inter-hemispheric connectivity, 
the predicted structural connectivity matrix had its largest 
eigenvalue at 0.87, which is in experimental agreements with 
EEG recordings of brain activity, and satisfies the theoret¬ 
ical and physiological conditions of dynamical stability and 
criticality in normal brain functioning. However, this solu¬ 
tion was not sparse and returns a full matrix, has negative 
entries, and shows a range of connection strengths different 
from the range in S, Eig. ih). 

5.3 Known structure recovery and prediction 
of inter-hemispheric links 

We now present the results of our algorithm applied to 
the brain data from [^. Eigures [^a-d] show the inferred 
matrices Ap, Xn, Xpt, and Xpn, respectively. Erom vi¬ 
sual observation, we see that both intra-hemispheric (strong 
main diagonal connectivity) and inter-hemispheric connec¬ 
tivity (strong off-diagonals) between the left and right hemi¬ 
spheres has been retrived. Eigures[^e-f] show the histogram 
of connection strengths for Ap, Apt, and Apn, respectively. 
The range of connection strengths inferred is very close to 
the distribution in the experimental DSI based structural 
connectivity matrix, Eig [^a) and (b). Our solutions pre¬ 
dict a higher number of both intra-hemispheric and inter- 
hemispheric connections as compared to the existing struc¬ 
tural matrix. The panel of plots on the right in Eig. 
shows the actual and predicted structural connectivities, 
along with plots of similarities and differences between ac¬ 
tual and predicted connections. The similarity plot estab¬ 
lishes that a high number of actually existing connections 
are recovered, as also confirmed by recall computation, and 
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Figure 5: Precision and recall with respect to 

varying number of dimensions preserved in spec¬ 
tral representation k for real brain data from [1Q| . 
Solid lines show precision, dotted lines show re¬ 
call, for Ap, Xpt with near zero entries removed and 
Apn = Ap n (An = 0). All have low precision, since 
original data has missing inter-hemispheric connec¬ 
tions and solutions identify this strongly. Recall is 
high for all solutions, with Apt showing best perfor¬ 
mance. 

the differences plot establishes that inter-hemispheric long 
range connections are retrieved. Note again, that we make 
no use of the experimental structural connectivity matrix to 
make our inferences, using only the functional connectivity 
data and applying our algorithm to it. 

Eurther, we also compare the results for precision and re¬ 
call in comparison with the experimental structural connec¬ 
tivity matrix. Eor all three matrices, Ap, Apt, and Apn, 
recall is high, with Apt performing the best at about 0.7 
- 0.8. It would appear that increasing the value of k fur¬ 
ther results in even better or near perfect recall. How¬ 
ever, we keep in mind the related important task of predict¬ 
ing the inter-hemispheric connectivity, that does not have a 
strong signature in the experimental structural connectivity 
data. Therefore, increasing the value of k results in the solu¬ 
tions moving closer and closer to the experimental structural 
matrix with decreasing intensities of entries for the inter- 
hemispheric connectivity. This observation also explains 
why precision is low for all three matrices Ap, Apt, and Apni 
precision is computed as the ratio of the correctly detected 
existing and mapped connections in the experimental struc¬ 
tural DSI data to the total number detected. Because the 
predicted interhemispheric connections retreived are compa¬ 
rable in both number and intensity to the intra-hemispheric 
ones, but do not exist in the existing and mapped structural 
DSI data, the precision goes down. However, we have shown 
in the previous section on the baseline synthetic model, that 
when no connections are missing, precision and recall return 
sufficiently high acceptable values. 

6 . RELATED WORK 

Brain network analysis is an active area of research across 
many disciplines including data mining and machine learn¬ 
ing, complex systems and of course neuroscience [^[^[7|[^ 
However, an optimization based perspective 














(c) (d) 



°0 0.05 0.1 0.15 0.2 °0 0.05 0.1 0.15 0.2 

Edge weights (structural predicted Edge weights (structural derived 'i^) 











0.08 


0.06 


0.04 


0.02 


0 

0.12 


Predicted connectivity: similarity with structural DSI data 


0.04 

0.03 Predicted connectivity: differences from structural DSI data 


0.22 

02 

0.18 

0.16 

0.14 

0.12 

0.1 

0.08 

0.06 

0.04 

0.02 

0.12 

0.11 

0.1 

0.09 

0.08 

0.07 

0.06 

0.05 

0.04 

0.03 

0.02 

0.12 


x-y coordinates x-z coordinates y-z coordinates 


Structural connectivity data from DSI 


Predicted connectivity total 


Figure 6: Results of applying algorithm on real brain data, (a) Xp, (b) Xn» (c) Xpt. (d) Xpn* (e), (f), 
(g) show histograms of connections strengths for (a), (c), and (d), respectively. Panel on the side shows 
comparison of actual and predicted connectivity. Columnwise: the plots represent 2D plots of nodes and 
edges with x-y, x-z, and y-z coordinates showing three different visualizations of the cortex. Row-wise: Top 
row shows the edges from the DSI structural connectivity matrix S, second row shows the edges from the 
solution matrix Xpt shown in Fig. [^c); note that the solution shows long range connections. Third row shows 
the edges from Xpt fl S', the part of S that is detected in Xpt. Last row shows the links that are predicted in 
Xpt but do not exist in S, note the long range connections. For ease of representation, solution plots have 
been thresholded at 0.02, so that too many very low weight edges are not visible. Our retrieved solutions 
have more connections predicted for both intra-hemispheric and inter-hemispheric connections. 


































has been mostly pursued in the data mining and machine 
learning community [^[^[^. These have mostly focused 
on reconstructions and inferences of functional connectivity 
from various types of functional connectivity data. 

A popular approach has been to represent neuroimaging 
or functional data as tensors. For example, fMRI data can 
be represented as a fourth order tensor consisting of the 3D 
spatial coordinates and time. Data mining tasks, like super¬ 
vised learning, have been extended for tensors for brain data 
analysis. In [^, data is assumed to be given in the form 
where Xi is the tensor representing an fMRI 
and yi is a binary label indicating whether the brain sam¬ 
ple is normal or suffering from a disease (like ADHD). Then 
an SVM optimization formulation is extended to learn the 
weight tensor of the separating hyperplane. 

In Davidson et. al. a tensor formulation is used to 
set up an optimization problem for both node and edge dis¬ 
covery. The motivation is that fMRI parcellation needs to 
be aggregated at the appropriate level which simultaneously 
corresponds to coherent functional regions while keeping the 
fMRI activity discriminative enough to capture local varia¬ 
tions. 

The body of work introduced in 2^ uses an optimiza¬ 
tion formulation to estimate a sparse inverse covariance ma¬ 
trix from the sample covariance matrix extracted from fMRI 
data. The insight for infering the inverse covariance matrix 
is the well known result that if data follows a multivariate 
normal distribution then the zero locations in the inverse 
covariance matrix corresponds to conditional independence 
between the variables, i.e., if 0(z, j) = 0 then i and j are con¬ 
ditionally independent and thus there does not exist an edge 
between them. Sparsity is enforced using the ii regulariza¬ 
tion on vec{0). However, the focus is still on the retrieval 
of relevant functional connectivity from data on functional 
connectivity. 

The inverse problem of inferring structural or anatomical 
connectivity is extremely topical in the computational neu- 

, but to the best of our knowl- 
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roscience community [7||19 
edge, has not been modeled with an optimization framework 
in the machine learning and data mining community. One 
of our chief contributions is that we show that there can be 
significant improvements in the quality of the predictions if 
we model the inverse problem using optimization. Our opti¬ 
mization formulation is also quite distinct in that we assume 
that the spectral representation of the fMRI adjacency ma¬ 
trix lives in a union of low-dimensional spaces where data 
points can be expressed as a sparse linear combination of 
other elements in the subspace. This has been referred to as 
the self-expressive property of data in , and we show that 
it is possible to derive an optimal sparse representation of 
structural/anatomical connectivity, starting from functional 
connectivity data using this property. 


7. CONCLUSIONS 

We presented a convex optimization formulation and ADMM 
algorithm for solving the inverse problem of deriving the 
sparse structural connectivity of the brain from its func¬ 
tional connectivity signature. To the best of our knowledge, 
this inverse structural/anatomical inference problem from 
fMRI data has not been modeled within an optimization 
framework before. 

Using only the functional connectivity data as input, we 
presented (a) an optimization problem that models con¬ 


straints based on known physiological observations, and (b) 
an ADMM algorithm for solving it. We showed that the 
algorithm not only successfully recovers the known struc¬ 
tural connectivity of the brain, but is also able to robustly 
predict the long range interhemispheric connections missed 
by DSI or DTI. Particularly, it addresses two principal gaps 
that remained unexplored in previous attempts: (a) our so¬ 
lution was sparse, and was very close to the sparsity that 
is actually observed in a DSI or DTI based experimental 
anatomical network, and (b) our predictions of the numeri¬ 
cal distributions of the weights or connections strengths also 
showed a good match with the DSI or DTI based experimen¬ 
tal anatomical network, even though we use only the fMRI 
data that has a very different distribution of weights. Fur¬ 
ther, we demonstrated these results on one of few available 
datasets that contain parallel maps of structural and func¬ 
tional connectivity of the human cortex parcellated into fine- 
scale 998 regions. Previous results have either been shown 
on the much coarser 66 node parcellation, and/or has not 
had a sparse form with close matches to the actual numeri¬ 
cal distributions of the experimentally observed connection 
weights. 

For future work, it will be possible to extend this formu¬ 
lation to include constraints coming in from the structural 
data and constraints that also respect the criticality and 
stability conditions that have been noted to be important in 
both and [^. For example, the problem could be for¬ 
mulated as a constrained matrix completion problem, where 
we write an optimization problem to use the fMRI data to 
introduce links into the DSI data in a constrained way, while 
also forcing the stability and criticality conditions on eigen¬ 
values as well as maintaining the required amount of spar¬ 
sity and the quantitative distribution of weights/connection 
strengths. The method could also be extended to model 
more realistic cases of asymmetric connectivity and/or com¬ 
binations of excitatory and inhibitory connectivity in the 
cortex, since in each of these cases, the basic optimization 
formulation could be changed to model different types of 
conditions. 

While the present paper focusses on brain networks, the 
inverse problem of deriving structure from function in gen¬ 
eral is relevant to many other types of dynamical networks. 
Infrastructure networks, such as roads or communication 
lines, have links that can support fixed flows. In this case, 
the structure of carrying capacity can be mapped exactly, 
but functional data on congestion and jams can be used to 
infer management or design strategies for the future, for al¬ 
tering the physical structure of the network, or altering the 
flow or capacities of links . As a second example, consider 
that while the exact structure of online social networks such 
as Facebook can be mapped, much of the structure may be 
dormant, as not all links are active all of the time. In such 
a case, deriving the effective structure of the network via 
studying time-stamped dynamics such as the frequency of 
interactions via specific links, can provide much deeper lev¬ 
els of information about the structure of the network. Thus, 
we believe that modeling this problem within an optimiza¬ 
tion framework as a general strategy can provide significant 
advantages and deep insights for several domains. 
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